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Abstract 


This paper iS conoerned with a computer program used tor Studies of the 
disturbed eoneS Around bodies in flowing plAsmaS, partiOv''ariy Spacecraft and 
their associated Sheaths and wakes. The program solves a coupled Poisson- 
VlaSov system of nonlineat’ partial -differential -integral equations to obtain distri- 
butions of electric potential and iOn and electron density about a finite-length 
cylinder in a plasma flow at arbitrary loh MaCh numbers. Using the Author's 
"insMe-out method"® which follows ion and electron trajectories backward to theit* 
origin at the body Surface or in the undisturbed plasma, together with a special 
iteration algorithm for self-consistency* the program takes into Account the parti- 
cle thermal motions with relatively few simplifying assumptions. The approach is 
applicable to a larger range Of parameters than Other available approaches. In 
Sample Calculations, bodies up to 100 Debye lengths in radius are treated, that iS, 
larger than any previously treated realistically. Applications are made to in-situ 
satellite experiments. 


1. mttOlMICTJON 


The problem bt theoretically calculating the structure of the disturbed piaswa 
(frequently referring to the wake and/or sheath) around a moving body m space in- 
volves the solution of a complicated system of coupled nonlinear partial differential/ 
Integral equations. I- The equations consist of the Vlasov (cOlllslonleSs Boltzmann) 
equations for the Ions and electrons, and the Poisson equation relating the electrlt 
field to the distributions of ions and electrons. The difficulty vS eesentially a 
numerical one because analytic Solutions are not possible (for cases of interest), 
and there IS no urtlqUe approach. In cases of stationary bodies. as iSell as 
moving bodies (theoretical references Cited by Parker 1). combinations of numeri- 
cal techniques (finite differences, iteration, quadratures, etc. > are required for 
treating Various parts Of the problem. For either Stationary Or moving bodies, 
the choices of techniques and their Use to achieve consistent solutions for any 
given set of physical parameters (defining body and plasma) have never been 
Obvious. Innovations are frequently required. The purpose of this paper is to 
present a technique suitable for a pillbox -Shaped body (v»ith emphasis Oh the 
^ake) 1 which appears tO be reasonably successful over a large range of the physi- 
d pardeufs. and » pr«e»t S««pla aolatlops IdcUldlne appUcaUand M 
spacecraft data. The pillbox problem Is illustrated schematically m Figur 1. 
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Various apiiroathes «hlOi have bson used tor this ln>e of problem are eum- 
marlaed bp Parber. a„ aubb calbulallobe. elmpu aaeumpflous a" Ido. 

ine custofnary ones aref 

(1) Collisions negligible (bul extensions of Parkei’-s Colllsional theory'* may be 
feasible iot the wake probleth). 

(2) Geomagnetic field negligible. 

(3) Simplified geometry (use 6f various types of Symmetries). 

4 Simplified surface reactions (usually, charged particles are neutralized). 

back /T' Phdtoelectron. 

scattered -electron, and secondary-electron emission are Includable).. 

ar. IludloT"”' «"“* 00-dbbtlvm.s 
Cl) steady state. 

dene^"? may be questioned (for example the neglect of time- 

dependent phenomena), but they may be at least partially relaxed by employing 
nown techniques to generalize the Calculations. In the interest of achieving 
r asonably ecohomical calculations Within the limits of available computers, the 
a ove assumptions in their usual form are adopted in the present work. 

The techniques and cortinuter program described by Parker* have been devei- 
oped to solve the coupled Poisson-Vlasov system of equations to obtain distribu- 

(with ° density, and potential, about three-dimensional bodies 

with axial symmetry about the direction of plasma flow). The method involves the 
use 0 a numerical grid or mesh of discrete points in space, with the potential and 
density distributions defined at these points. The Poisson and VlaSov equations are 

cretTltm ‘ gHd PdtntS. A sample of such a dis- 

rlZ T\ 'ir' T"" 

cvlinae Associated with each point is a volume, ih the form of a 

cylinder for points along the axis, ahd in the form of a torus of rectangular cross- 
section for all other points, ccmnguiar cross- 

Jw ^ Sbftbce Is shown b, s hssvy oumne in tbs interior of tbe 

Lnin I .h Md columns of jrtd 

po in ts s s Shown. HC^. grid point, are un«,«all, Spaesd. so tbs. . higher 

density of ^ ^ ^ 

The nolemM 7’'".“*°*’ ° “*"< "ttlelMUy 

probe ) while the remainder IS at another polentlsl. Ifhe surf.ce cun con.1.1 of 
^rilons with srhltrurlly uesl,pied eonducllvlly Sbd emiuelbb churucllell ft. 
actual number of grid polnte used Was of the order of bnudrede. rather than tens 
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Figure 2. Discretizatloti in r'-t Space 


aS Illustrated. ) The shaded areas Surrounding giMd points are the croSs- sections 
Of toroidal volumes aS in Figure 2. At the outer boundary of the grid, one must 
represent numerically the boundary condition at infinity, namely, SUCh that the 
potential vanish and the velOOity distribution be the unperturbed one. This boundary 
must be sufficiently far out to represent the outer condition accurately. It turns 
out to be more efficient to Use a "floating" rather than a "fixed" condition On the 
potential (Section 3). 

While the present problem is axially symmetric, it ean be generalized to three 
dimensions aS follows. The grid in Figure 3 consists of points in r-z Space, and 
the assooiated volumes are tori. The generalization would Oohsist Of inoluding the 
aSlmuthal variation by adding an aSimuthal angle 6 to the coordinate system. The 
disoretization In d would consist of having a number of azimuthal planes in r and z. 
each labeled by a given value Of d . Thus, for example, the r-z plane of P*igure 3 
would be characterized by a given value of d. The volumes associated with the grid 
points would then be pie- shaped. 

In the next section ^Section 2) the "InSide-out" method for evaluating particle 
fluxes afad densities (solving the "Vlasov problem'*), developed by the author in 
lP64, will be discUSsed< with reference to the grid of ^'igurS 3. In this method 
the ton and electron trajectories are followed backward in time, from the point Ih 
space at which it is desired to khOW the Velocity distrihutioni to the source of the 


^8 North boundary 


probe 


sotellite 

eurface 


South boundary 

Flgurfe 3. Dlfferfenc6 Equation Grid 


particles in the undisturbed plaSma or at the eurface, where the dietributlotis are 
known. Figure 4, referred to again in the nekt Section, illustrates Schematically 
hoW a trajectory is traced backward from any point P and is found either to reach 
the body surface or "infinity*' at the boundary of the grid. (The point P is usually 
but not necessarily one of the Spatial grid points: it can also be a surface point. ) 


The "delta-factor" is a cutoff function, and i, J, and k are indices associated 
with one of the trajectories used to evaluate density or flux as discussed in the 
next section. The case illustrated is for contributions from the ambient plasms; 
for contributions from the surface, the values of 6 (zero and unity) are Interchanged 







Evaluation of for (ij.lc)-th trajectory by following (reversible) 
trajectories backward in time. 

Figure 4. Basis of the Inside -Out Method 


Figure 5 illustrates the four possible types of trajectories which can contribute 
to the particle density at a point. These are Types 1. 2, and 4, so-named by 
Parker^ and defined as follows: 

Type 1 

One-way trajectories, going from infinity to the surface, or from the surface 
to infinity. 

Type 2 ^ 

rivo-way trajectories from Infinity, which cOfne in, pass through a position at 

mlnlinum distance frotfi the body surface, and go out again. 

type S 

IVo-way trajectories from the body surface, which go out, pass through a 
position at maximum distance frOm the body surface, and come in again. 

Tj^e 4 

closed or nearly- closed trajectories which orbit about the body indefinitely, 
these can only be populated and depopulated by collisions, which are neglected in 
the present work. An analysis of the effects of coUlsions on type-4 trajectories 
has been performed by Parker. ^ 

it should be noted that contributions to surface flukes can be comprised only of 
type- 1. and type-S trajectories* while all four types contribute to sp>ce charge. 
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In Section 3 the "Poisson Pfoblem" is discussed, whfere the electric field 
(potential distribution) is computed with the ion and electron densities considered 
known.. On the other hand, the "VlaSOV Problem" (Section 2) involves computation 
of the ion and electron densities with the field considered known. Hence, since 
neither the field nor the particle densities are known initially, the Poisson and 
Vlasov problems must be Solved Simultaneously. 

An iteration method may be used for computing Self ‘■consistent charged* 
particle and potential distributions, this is herein referred to as the "Poisson* 
viasov iteration. " two principal options are employed for this procedure in the 
present program. In one of the options, the "charge>density" option* the space 
charge is initially and arbitrarily assumed to be zero. For this case, one obtains 
the Laplace (Space-chafge-leSS) electric field from the t^oiSSon problem, this is 
the "zero-order" potential distribution, which becomes input to the Vlasov problem, 
the resulting solution of the Vlasov problem yields the ion and electron densities 
at the grid points, which are combined to mal-j "zero-order" charge densities, 
these become input to the hSkt PoiSsOn problem, which then yields the "first-order" 
potentials, and so on. in this procedure one usually "ml*es" successive charge- 
density iterates to improve stability; otherwise, the process can "blow up. " One 
can also mix potential iterates rather than densities if desired, the dependence 
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of the atabllity and cotiVergcnhe 6t die hboi^e procedure On the mticliig parhnletei^ 
have been studied analytically by Parker^ and Parker and SulllVan. ^ (No other 
analysis of this type has been published to the authoi:‘'S knowledge. ) This (charge- 
density) option is most effective ivhen the spatial region of interest Is not too many 
Debye lengths across. The analysis ShoiivS that one can (probably always) choose a 
mixing parameter Sufficiently Small to ensure convergence, but at the expense of 
additional iterations. 

In the other Option, the "ion- density option, " the ion density distribution alone 
is assumed initially. Initial guesses which can be employed include (1) aero ion 
density everywhere, (2) unit ion density (the ambient value) everywhere, and (3) the 
neutral ion density which obtains When there are no forces. Whichever choice iS 
made for the initial guess is designated the "zero-order" ion density. Now if one 
can assume the electron density to be given by the Boltzmann factor exp(0), thus 
avoiding trajectory calculations for the electrons and affording computer economy, 
the Poisson equation may be solved, holding the ion densities fixed, but regarding 
both the potentials and the electron densities at the grid points aS unknowns. This 
is a nonlinear problem, which is solvable by a modification Of the relaxation pro- 
cedure used for the "charge-density" option. The new procedure is art important 
advance Since the iteration ia not as sensitive (tending to blow U0) fOr small Debye 
numbers aS in the charge density option. Thus, very large bodies (in multiples 
of the Debye length) can be treated. This has been the method used to obtain the . 
large-body results shown below. Similar ideas have been used by Call and 
Fournier, but these workers have not treated large bodies. 

The assumption that the electron density is given by the Boltzmann factor 
becomes invalid when the body surface potential is near zero, or when there is a 
potential barrier or "well" in the wake Such that the wake potentials are more nega- 
tive than the Surface potential (causing electrons to be attracted to the surface 
rather than repelled from it). In this case it IS still possible to use the ion-dCnSity 
option, with its large-body capability, provided that, within each "major" iteration 
cycle a "minor" iteration is carried out With the ioh densities held fixed such that 
the electron densities are computed realistically by trajectory calculations, at 
least for points near the Surface. 

This latter technique is as yet in an experimental stage, but it seems promis- 
ing in that it may produce Solutions with reasonable costs for large-body problems; 
in such problems, the conventional PoiSson- Vlasov iteration based On the charge- 
density option becomes expensive. ^ A disadvantage of the ion- density option, 
however, is that its convergence properties are not understood; therefore, its costs 
are difficult to predict. This Is in contrast to the case Of the charge- density option 
where an analysis is available. ^ 
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Before cbnstdeiilhg fUirther-detatlb, we mbke here soi»e general remarka cbn^* 
cernlhg the methbd. FolldWing thia» the principal results will be sutnihariaed. 

Briefly, the present approach^ differs fpom these of Call® and Martin^* by 
Including both the loh and the electron thermal motions, whereas Call and Martin 
represent the distribution of ions by a cold beam and use an "outslde*lh" method. ^ 
The approach differs from that of Taylor^^ in that (1) it is applied to three-dimen- 
sional bodies wher>>as Taylor treats an infinitely- long "thick Strip" of rectangular 
cross- Section, and (^) die Poisson and VlasOV Calculations are cycled Until Self'- 
conslstenCy is achle^^ed, whereas Taylor's calculation is not Self- consistent because 
it is terminated after the first cycle. The approach differs from that of GrabowSki 
and Fischer because they (1) assume that quasineutrality holds everywhere (an 
invalid assumption in the very near Wake - See below and Section 5), and (2) apply 
their method to an infinitely-long cylinder. Differences with other methods are 
outlined in Parker. ^ The most similar calculation previously done Was for an 
infinitely- long cylinder by Fournier, using the inside-out method. The present 

author has used the method for two-electrode rocket-borne and laboratory probe 
4 ■ A 

systemj, for the problem of a small probe in the sheath of a large electrode, and 

most recently for the problem. Of die pillbox- Shaped Spacecraft. ^ The inSide-out 

method Was also uS6d by Parker and Whipple for the theory of a satellite BuSh- 

moUnted prObe. 

Two major advances are represented by the present program, as opposed to 
previous- approaches, particularly with regard to Wakes of three-dimensional 
bodies: 

(1) Thermal motions of ions as Well as of electrons are treated realistically 
by following their trajectories in the electric field. (The ion and electron tempera- 
tures Can be different. ) 

(2) The technique for achieving self-consiStenCy is promising for large bodies 
many orders-of-maghitude larger than the Debye length (the iShuttle-Orbiter or the 
moon, for example). 

Solutions may be obtained With reasonable amounts of computer time by 
judicious choices of grid points and other numerical parameters. The method can 
be extended to include an arbitrarily- shaped body (presently a body of revolution), 
electron emission from the Surface, and differential Charging when the surface 
consists of Sections with different conductivity. 

In Sections 2 and 3 some details of the techniques for the fluk and density 
calculation ("Vlasov problem"), and for the PoisSOn calculation, are treated. 

Sample calcUlational results are presented in Sections 4 and 5. 


t 


1.1 Summary of Principal Results 

The prtocipaljresults aire as. follows, tn-situ experimehts associated wiUi. the 
AHel 1 and Explorer 3 1 satellites Si-e modeled by a pillbox geometry. The Artel 1- 
ei^erlment obse!r«es distinct Wake structures ai^Sociated With die main body and a 
Small external Ion probe. Transverse profiles Of electron current are measured 
at S main-body radii downstream. The two wake structures are similar in ttiat 
they both show a below -ambient Central core or peak within a depleted-region of the 
order of the width of the main body. The theoretical results for the assumed values 
of the parameters associated with the experiment show no well-definsd central core. 
They furtiier show that all structure dies away beyond 6 or 8 radii downstream. 

A pronounced electrical focusing of ions in the wake of the highly- negative ion probe 
is predicted by the calculation, but this disturbance is Confined to tixe relatively 
near wake and does not persist downstream. The filling of the wake in bodi cases 
by tiie plasma suggests a-fluid-llke bulk motion of the plasma. As one moves down- 
stream with the plasma, the motion is at first radially inward, this is followed by 
a pile-up and a single "boUnCe" after which the motion. iS-oUtwacd. — Simultaneously, 
the disturbance becomes weak and dies away. 

Two Explorer 31 cases are computed for different values of the ion Mach 
number (the Other parameters remaining roughly comparable), and in both cases 
the body is-several kT negative. In the Case of the larger Mach4tUm:-er (3. 4), the 
ion density in the near wake Is beloW the Corresponding electron density, emd both 
are significantly below ambiet^ This is consistent with the traditional picture of 
wake structures with ion MaCh numbers significantly above unity. In ^e case of 
^e lower Mach number (l. l) the ion density in the near wake is above the corres- 
ponding electron density, and moreover the ion density is roughly ambient. This 
latter may Seem unexpected, but is understandable on the basis Of Langmuir-probe 
sheath theory: In me sheam of a slowly-moving negative probe me ion density 
should predominate over the electron density. This latter result is new in wake 
meory, and arises because Ipw-Mach^'number wakes with.space charge have not 
been previously rigorously Computed, 

For me wake of a large body (lOO Debye lengms in radius) in me form of a disk, 
me results show mat quasineutrality is valid outside of a cone-shaped region in me 
very near wake, and is invalid wimin mis region. Other features of the large-body 
wake structure include (1) a potential well in me near wake, and (2) a central core 
of approximately ambient density of both iOns and electrons. This latter feature 
seems similar to mat observed in me Ariel 1 experiment. 
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2. Till': iVi^iii^'dlit-MK'uioi) 


Tliere is vtidts than otte a|)|)rohch t6 th6 problem ot chlcul&ting sheath aiid W&ke 
structures, these approaches have In cOtmhdh the following elements, the quanti- 
ties to be computed include (1) the potential distribution, and (2) the iOn and electron 
density distribittions. One may also include the associated Surface current densi- 
ties. The equations to be solved simultaneously are (1) the Viasov equation for ions, 
(2) the VlaSoV equation fOr electrons, and (3) the PoiSSOn equation. The solutions 
of the Vlasov equations (velocity- distribution functions) are USed to compute number 
densities (and Surface Current densities). The number- density distributions become 
input to the (right-hand side of the) Poisson equation Which yields the potential 
distribution. Finally, an iterative procedure is used for self- Consistency, wherein 
the density and potential distributions are successively Cycled Until satisfactory 
convergence has beeU achieved. 

The steady- State Vlasov, equations for ions and electrons state that the velocity 
distribution, functions remain constant along particle trajectories. With the electric 
field assumed given (numerically in terns Of a Spatial grid about the body). Solving 
the Vlasov equations means formally that one determines, from the shapes Of the 
trajectories, the ion and electron velocity distributions at the grid points. The tra- 
jectories relate local velocities at a given grid point to those at infinity or the sur- 
face. Through these relationships, the iOn or electron number density at the point 
may be evaluated by a velocity-integral over the local velocity distribution. Sim- 
ilarly, the current density may be evaluated at desired locations (usually the body 
surface). 

It is convenient to classify various theoretical approaches on the basis of hoW 
they treat the trajectory part of the Vlasov problem. Ail "inside-Out" method fol- 
lows the trajectories backward in time to their source, while an "outside-in" method 
follows the trajectories forward, in the direction of physical motion of the particles. 
(Di an oUtside-in method, the velocity-distribution function is not calculated; rather, 
the density is evaluated directly. ) There are in addition other (less realistic) 
methods involving approximations where trajectories are not followed at all. The 
three types of approaches are discussed In Parker. ^ There exists as yet no 
systematic comparison of the results of the various approaches With one another. 

^’or the purposes of discussing the thslde-out method. We define here the para- 
meters Of interest: 

ftasma Parameters 

n^ = unperturbed number density at litfinity 

^1' ®l^ctren temperatures (= T for SqUai ion and electron 

temperatures) 

m^ = ion mass (electron mass not required) 


I 


» electron Debye length 
Body Par^eters 

B. - cheracterUtlc diitiebslon 
o 

= relative velocity of body and plasma 
= body potential 

^0 “ “ dimensiOMlese body potential 

M = V Vm./2kT. = Ion. Mach number (electron Mach number assumed 

® ^ ' negligible) 

Xjj = = Debye number 

Henceforth, all lengths are to be considered normalized by B^. Thus, X^^ will 
denote the dimensionless Debye number. We also consider here the Case of equal 
temperatures. Potentials are normalized by kT/e, so that (^(r) denotes the dimen- 
sionless potential at the Spatial point r. Number densities are normalized by n^. 
so that n( r ) denotes the dimensionless density at r . In the calculations involving 
integrations Over velocities, v will denote a velocity normalized by the value of 
associated with the particles 6f interest. Similarly, E will denote total 
energy normalized by kT. Velocity- distribution functions (denoted by f) will be 
normalized by n^. t'or a given body geometry, there are three dimenSionleSS 
physical parameters of interest, namely, Xj^, and M. (Eor unequal tempera- 
tures. the temperature ratio T^/T^ represents an additional parameter. ) 

Consider a single species Of (charged) particle, that is, ionS or electrons. 

The electric field is assumed to be known. In order to compute the number density 
n(r ) at the point r, one must evaluate the triple integral over velocity space; 

n( r ) = ^ **^y 

Where t(r,^) is the distribution function which satisfies die Boltzmanh equation 
for the given species Of particle, r is the radius vector of the Space point Of interest, 
and V is the local Velocity Of a particle at r. The velocity-volume element is 
written as if cartesian coordinates were being used, but the product dv^^dv^dv^ is 
intended to symbolize an arbitrary coordinate system. Similarly, in order to 
compute the collected flux at points on the Surface of a body, one must evaluate at 
each point a triple integral over velocity space Of the form 



' 2kT/m 



Whel-ejTjj U the component of the pafticle veloelty normal to the sorface at the 
pdiht r. 

The problem la thus to evaluate f. Since the problems Of interest are assumed 
to be colllslOnless and constant in time, the distribution function f satisfies the 
Steady-state VlaSoV (or ColliSionless Boltzmahn) eejuatiort, namely, 

v-=Vf+a. = 0 ^3, 

afhere a is the vector acceleration of a particle passing with Velocity v through 
th^point p. The gradient operators Vand operate on the Components of rand 
of V, respectively. Equation (3) states that f is constant along a particle orbit, 
which is characterized by the constants of ihe-motion. In a general electrostatic 
field (here assumed given) whose sources are Volume and surface charges, the 
total energy E is conserved, where the dimensionless E is defined4>y 

E-’.^ + «(?) (4) 

and Hr) iS the dimensionless potential energy of the particle at r. 

Witii *(r) a known function of r, one may evaluate the integrals in E<js. (1) 
and (12) by following orbits backward in time with trajectory calculations to a point 
where f is known. For example, in tiie case of a body immersed in a plasma, f 
is assumed to be known at infinity (where <l> vanishes), and is assumed to have at 
infinity a prescribed energy distribution, such as a Maxwellian with drift, or a 
more general distribution. Also, f is assumed to be known on the surfaces of the 
Spacecraft. If a surface emits particles, its distribution function must be pres- 
cribed. If the Surface absorbs without reemitting Charged particles, the distribution 
function (of emitted particles) is prescribed to be zero. Thus, f is discontinuous 
in velocity space. That is, the physlcally-possible velocity Space (at the point r) 
is divided Into three domains, namely, the domain of orbits which have COme to r 
from infinity, the domain of orbits which have come to r from the spacecraft Sur^ 
faces, and the domain of trapped orbits (assumed to be unoccupied). Ihe shapes of 
the boundaries between the domains depend, of course, on the geometry and the 
potential function ♦, and it is the heart of the problem (l) to determine the bounda- 
ries of the domains of orbits, and (2) to evaluate the Integrals £qs. (1) and C2) over 
those domains of velocity Space. 

In practice, orte need not in general determine explicitly the boundary of a 
domain in velocity spaOe. Rather, one may follow a large number of orbits back- 
ward in time ^computationally), and evaluate the moment integrals^ Eqs. (l) and 
(2), automatically from the results of the orbit- following, it may, however, under 
some circumstances be more accurate and efficient to determine this boundary, tb 
db 80 would complioaie the cbmputer pr6grcuhtning» 
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M ‘"••"““to" no.i«g «» »-di«cu„„ ,ttk 

M, the dlMenstonleae eeloolty-aielHhUtieh ftlhetta, >t totinlty me, be «-MtteO! 


1 

,3/2 


-(v* + M* . 3Mv^J 


1 -W + v^ + M^ - 2Mv^ » 
i7? e *<* 


.3/2 


( 5 ) 


(velocities in Uitits ot 2kt’/m , 

= axial Coiia|)onent of velocity) 


I 

>«h6re V ' 


^ ♦ may be identified with the total energy E. and w ^ witht/i" 


ribution may also be used to represent particles emitted from the surface, tlie 


■ ■///■ 


6t 


‘’''“5S|Ayhj'lJkVuk 


] 


( 6 ) 


e3- 


•here d j le a ehort-hahd notation for the eleineht d» d, d, and » te a edtoH 

tor etep) tonetlon. e,ual to onlt, Or aero, depend.,* on (1^ .Lther the trajeetor, 

found tft . ** '-'■wx jf 


' • o A*/ wax^^fcucx- Liie iraiecrorv 

.0 teond to oon.e fron. idhttty or the bod, enrf.o.. and (2) .hether n repreeent. 

^ Tl ^ “”•■^'1 »r pertlelee irom Inflnlt, or Iron, the nnrtaee. In the eon.. 
^ three Indtoea refer to dUcrete *aluee of three cooponents of veloott,. .here 
the ranee a« chosen In aooord .Ith a quadrature seheme (Gaussian), and the 
coefficients .re pnoportlon.1 to the associated .eights and other factor.. Bach 
to™ fl “ ‘“'•Wdual trajeetor,. a .hnllar Sum is obtained 


In order to erauate Eq. (g) (or the denstfy. or the corresponding equation for 

" :ier rrr " 


& « polSr aiigle With resfiect 16 e-axis 

0 ‘ aeiinupsi angle with respect to the plsne contSihing the e-axis and the 


point r . 


The angles a and 0 which define the orientation of the VeWdty-vectbr v (whlla./f 
defines its magnitude) Sre illustrated in Figure 6. ^ iwntie^E 

It is shown by Parker^ that the Integrals for both the density end flux eSn be 
transformed to the following form suitable for QaUssian quadratures: 
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Figure 6. Angle Variables iit Velbclty SpSce 


I 




6 (&, b, c) • da db dc 


(7) 


Hbre, the ranges of o,P. Snd 6 haVO been trsneformed to the Interval (-1, 1) through 
the use of new variables a.b, and c. In terms of these latter variables, o.(3. and 

E are given by 

er(a) * cOS'^a fOr density 
ff(a) = sin" ^ V ^4-^ for flux 



A' 


4 


‘ 

1^1 

1 .1 1 1 

i 1 < 1 


^(b> J (1 + b) 


■ fH-+ 


(9) 

( 10 ) 


Where the range of a is from 0 to » for density and from 0 to »/2 for flux; the range 
of 0 is from 0 tb «f; and the range of E is from ♦g to infinity, where denotes the 
potential of the Sdurce (Infinity dr the Surface). 

The-deflnltlon of T in the integrand of Eq. (7) is as follows: 


T(a.b, e) 


•U(a,b, c).. 


(1 - c)2 


for density 

C(-. 

T 


E(e) - ^ 


( 11 ) 


for flux 


where U(a,b, e) is given by 

U(a.b.c)= - «g + E(c) + M2 + 2M>^(^. coso^(a.b.c) 


( 12 ) 


'*^i*** ®oo denoting the value of the polar angle of the velocity- vector at infinity, 
whfoh depends on the local values of & , 0, and E (through a, b. and c). Hie product 
yE . coso^ in Eq, (12) is identical to the quantity v^^^ in Eq. (5). 

Now it is convenient for flexibUity to divide the a-doraain into equal sub- 
intervals, the b-domain into Mjj equal subintervals, and the c-domain into 
equal Subthtervals, and then to use Gaussian quadratures of order 2 in each 
subinterval. This leads to a sum of the following form : ^ 


2 
a 


M. 


M 


I-A..SS 


M^M- M 
a D e 


v V y; V V y 

^ JLj T(ajb;e')- 6(a!bJc') 


*^a=^ *^b=l '^b*l Je“^ 


where 


(13) 


- 4. / (-1)'^° 


+ 2K - 1 - M 
a 


.) 






*b\v^ 

= -l(±iL 
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(14) 


+ 2fc - 1 . M 
e e 
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Finally, the form ot fiq. (fi) may be obtained fram Eqs. Ii3) and (14) by 
Writing 


I 1 J R 

^ " ^IVe ? T(aj.b^.c^)- fi(aj.bj.c^) 

where 

a, = a* with i = J + 2 (K - 1) 

^ & Q 

bj = b' with j = Jj, + 2 - 1) 

c. = c' with k + 2 (K - 1) 

• e e 

I *= 2M J = 2M^ K = 2M^ 

Si D 6 


(15) 


(16) 


The SUhi of Eq. (13) consists of 8 terms, each bf which represents a 

trajectory, followed backward from the point of Interest. The cutoff function 
6(a^, bj, Cj^) is aero or unity, depending upon the trajectory defined by the indices 
i, j, and k characterizing the Initial velocity components. The caSe illustrated in 
Figure 4 is tbt contributions from the ambient plasma; for contributions from the 
surface, is aero Or Unity according as the trajectory reaches infinity or the 
body surface, respectively. 

The computed results to be presented later are baSed on the assumption of 
no surface emission. 

The method of computation of orbits involves integration of the equations of 
motion, with the forces given by the components of the gradients of potential. These 
components are obtained by interpolation between Values of potential defined at the 
points of the grid. Say of Figure 3, as described in Parker. ^ The criterion for 
"escape" or "absorption" of an orbit (that is, evaluation of 6) depends on the 
geometry of the problem and of the grid. The equations of motion are integrated 
step*by*step until the orbit either passes out of the outer boundary of the grid 
("escapes") or returns to the spacecraft surface (is "absorbed"*. The orbit compu- 
tation time- Step is not of physical importance in these time- independent problems 
where only the shape of the orbit matters. The time-step is kept as large as 
possible consistent with maintaining the energy loss or gain within desired limits. 
The method of integrating the equations of motion, the interpolation method to find 
the forces, and the control of step size, are discussed Irt Parker. * 

An important consideration is the accuracy of the quadrature- sum. Naturally, 
the accuracy is related to the number of terms used, that is, the number of orbits 
where each term corresponds to a unique orbit, in a test of the energy quadrature 
alone, and with M = O, the unperturbed value of density (unity) was computed for 
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, . A# wi » I 2 4 8 16, and 32. The cofresjidndlng nuirtertfcal enfore were 
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aCrtlrae, Inlthln a ,a» pafM«t> to. porpa... M computing 
ttaxwelllan Olatrlbutlon flthout OHft (or, for alertfon.1, for Urge ' 

(M) the accuracy of the above unmbdlfted quadrature is diminished. Modlficatlo 
Lr Improving th. accuracy at largo M by .ultably uclghtlhg the intogrand in the 
domains of importance are given In Parker. 


3. 


THE POISSON PROOLEH: POISSON DIFKERMCE EPUATIONS 


In the present problem the electrostatic field is axially symmetric and is 
defined on a meSh of spatial grid points, such that at any point (Incbiding gr 
oolnts) the potential and electric field can be obtained by interpolation. 

Assume that the Space charge density is known at the grid points. Consider a 
grou^ rlrlor grid U». tormtn, a portion of tbo overaU grtd a. abo.^ 

FlguL 7 in this figure, the vertical and horizontal directions are the a and r 

«(b feanectlveiv Where a and r denote the cyllndrlOal axial ahd cyUndrical 

rldtal'tooUlnatca. raap’octtvaly. Utrcc bortcontal 7‘“d“ 

* z and Zt^,. and three vertical grid lines, of constant rvalues rj.j. r^. 

are Shown In the figure. (Note that the index (1) of z increases as ^ > 
•i^e set of grid lines Intersect at 9 grid pOlhtS. or nodes, as Shown. Each poln 
he considered to be associated with a volume of Space, and to ^ve a group^ 
four neighboring points which "interact" with It. Thus, consider the central po nt 
of the gLp. labeled C in the figure, which may be identified with one of the grid 
Lint, to ftgurc 3. Assocteted »tth thl. point 1. c Volume of rcvolutton (. torn.) 
nboiC Cfoe.-..ctlon 1. rectcngular and 1. .hown by the rectangular .haded area 
iLLuX Lint C. Ute .haded area U defined by connect.^ the mtdpolnto of 
rTunou^dln, me.h rectangle.. Let V 

the neighboring point, (rtiove, belo«, to the right of. and to toe left of c) 
lUmled «. S. E and W (north, .onto, ea.t and *e.t, reipeettvely). 

Let the Poisson equation be written in dimensionless form aS 


« .p B (ttg - 


(17) 


M.i.gato.3 « n X ^ and d denote the dimensionless electron density, ion density. 
X :Shl.Vcto:.Lic potential and apace-charge dennity. feepeCweiy. and 
all lengths are in units of the body radius. 
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Figure t. Group of Interior Grid Points in r-i Grid 


The grid lines may be 
intervals are nonuniform, 
obtained by integrating Eq. 


considered to be arbitrarily chosen so that the mesh 
In this case the PoisSOn difference equations may be 
(11) over the volume t of the torus associated with 


Point C: 



^dt 



(l6) 


» Shown since t to omdl in principle, end pj, '» “■» " “* ^ 

Jlvergence theorem, , the left-hand side becomes 



( 19 ) 



where £ denoted the surface of the torus; dd/dn Is the component of vd in the 
outward normal direction at the surface; Ag. Ag, and denote the areas Of 
the north, south, east, and west surfaces, respectively; and the quantities 
(8^/8n)j|, g B vv Values Of B^/Bn taken to be constant on the corresponding 

surfaces. 

<8^/8n)jj g B W approximated by difference quotients, namely. 



<Vl-*l^ 



(dw - d) 

Of j ^ . ..■ 

<h - *i+l’ 



(»E • »> 

Vl ■ 



* ♦) 


( 20 ) 


where denotes the potential at Point C and dg. ^b* *1®*'®*® ^® neighboring 

potentials. If Point C iS an interior grid point, the areas Aj^. Ag. Aj,. and A^ 
are given by 


f ^ *‘j-l^ ^*i-l ■ *1+1^ 

and the volume t id given by 


T 



*1+1^ • 


( 21 ) 


( 22 ) 


thus we obtain the difference equation in the form 
where 

c = Cjj + Cg + Cg + 


(23) 


(24) 


and 
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" ^1+1^ 



, ^w 

r 


t 

I 


(25) 


This shows how to fom. the diffei'ence equations used fbr the f*blSsoA fiirableins of 
this papei*. Equation (24) ho’ds only tor an "intel^loi*" point of the grid, that Is. a 
point Surrounded by neighbors On all four sides. If Point C has a known neighbof* 
ing potential (fon example, if Point C is adjacent to the spacecraft surface), then 
the corresponding textn on the left*hand side of Eq, (23) Is transferred to the 
rlght'hand Side as a known quantity. 

The boundary conditions for the potentials In the Poisson problem are as 
follows. At points representing the body surface, the normalised potentials are 
fixed at the chosen values. At the external (boundary) points of the grid, where 
"infinity" is represented on the computer, a "floating" condition is optionally used, 
namely, a linear relation between ^ and 8^/dn, the normal component of The 
exact relation of $ to d^/8n is not Important when the external boundary of the grid 
Is sufficiently far away. (For the calculations to be reported, the assumed rela- 
tion was the same aS for a Coulomb potential. ) In any case, either the fixed 
condition ^ = 0 or the floating condition will give the same results, provided the 
grid boundary Is moved sufficiently far out. The effect of various types of boundary 
conditions representing "infinity" have been studied by Taylor and by Parker 
and Sullivan. ® in general, the floating condition appears to be computationally 
more efficient than the fixed one. Of course, the floating condition becomes ideal 
when the true relation between 0 and 8(fr/8n is used, but this requires that the 
asymptotic form of the solution be known i,i advance. For example, see Parker and 
Whipple. The boundary conditions at the outer grid surfaces can be combina- 
tions of fixed and floating conditions. 

Consider a Point C on the outer boundary of the grid where a floating boundary 
condition is chosen. If the potential is assumed to satisfy the linear law 


M = -tt4> 


8n 8a 

on the z-boundary (North or ^uth), and 


(26) 


B± 

dn 


(27) 
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dh the i^xboundary (Eadt only; Oh the Westy, then the correflpohdlhg "noighboh 
tetm" on the left-hand Olde Of Eq. (23) Vhnlshes^ahd the corhespondihg "neighbor 
coeflictent" on the right-hand aide of Eq.. (24) la replaoed by oA or 0A, where A ta- 
the appropriate area. The quantitied a and 0 depend Oh the position and oh the 
aasutned model for the variation Of the potential at large distances. 

Once the coefficients of all of t’le equations (corresponding to the grid pointa 
where the potentiale are unknown) are computed, the ayatem of linear equations of 
the form of Eq. (23) may be solved by iteration. Point- successive Over-relaxation 
is a well-known process and has been found to be effective in the present problem. 
For the relaxation process, one rearranges the equations. So that the "diagonal" 
term is alone on the left>>hand side, while all the other terms are on the right-hand 
side with the known charge -density term. Thus, Eq., .(2gi become# 

Cd = Cjj#jj_^Cgdg + + Pj,T . (28) 

First, an initial gueSs is made for the values Of all the potentials. Then new 
values are obtained from the left-hand sides of all Of the equations (28), using pre- 
vious values on the right-hand sides. One "sweeps" through the equations Success- 
ively. replacing the potentials oh the right-hand sides wlth^lpdated-values as they 
become available from preceding equations. This procedure is Usually Stable ahd 
leads to convergence. "OVer-relaxatlon" is the process of mixing Successive 
potential iterates in Such a way as to enhance the rate of convergence. ^ 

When the potential distribution is such that the electron density n^ is apprOxl- 
mablc by the Boltzmann factor exp(^), the relaxation equations (28) can include the — 
electron density as an unknown function on the left-hand side. The equations are 
then nonlinear; they may be solved for d by a Newtonian process, with the ion 
density n^ considered fixed. This procedure is promising for large-body problems.^ 


4. SAMPLE RESULTS APPLIED IN IN-SlTU DATA 

The calculations reported here refer to two in-situ experiinents, Ariel l^^ 
and Explorer 3l^^, where data are available. These results are preliminary in 
that they are intended as an illustration Of the capability of the program rather 
than as a systematic study. Geometrically, the body is assumed to be a pillbox, a 
cylinder Of height equal to its diameter. The three dimenSiohleSS physical para- 
meter# defining the problem are M, and defined (earlier) by: 
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M- s_ ,«ivJ‘/2kT 

whfere T la the pUfiiKia temperature (aasumed to be the aame for tons and electrons). 
^ IS the body potential (for a conducting nohemitllng body), la the velocity of 
the plasma flaw relative to the body and parallel to its axis, m la the ion mass, 

«o is the body radius, and Xj^ is the dimensional Debye length. 

The numerical parameters for the calculations to be described include of the 
order of 100 grid poims. distributed mostly in the wake region, and of the order 
of SOO trajectories per grid point, distributed among the two angles and the energy. 

4.1 Ariel.! SalellU^-- 

Figure 8 is a schematic drawing of the Ariel 1 satellite, showing the location 
of electron and ioruprobes. after Henderson and Samir. The- boom -mounted 

Sipis Silf 


tea 







Orobe rtieasures electpoft-currenta at a distance 5 from the Centet* of the satel- 
lite (ttialn-body) Which has a radios = 11. 5 In or 29 cm. The loh probe mounted 
near the surface and on the spin axis Is a small sphere 6. 6 times Smaller than the 
main body. The satellite oeloclt; is such that the Ion Mach number is about 4. The 
satellite potential Is about 4 kT (1 volt) negative with respect to space. The satel- 
lite radius is equal to about 10 Debye lengths. Due to the satellite motion, spin, 
and .orientation, the boom probe sweeps through the wake during each spin revolu- 
tion. In successive revolutions, it sweeps through at different_angles and samples 
the structure of different parts of the Wake. 

Figure 9 shows normalized electron current data taken, from the paper by 
Henderson and Samir (their figure 4). In. particular, the data at 0 = 84® aabeled 
"MAIN") samples the wake structure associated with the main body,. While the data 
at 0 = 60® (labeled "I. P. ") samples the wake structure of the ion probe. We will 
Consider- Separately the main body and ion probe in-the-foUowing comparisons 
between the data and tluaoretlcal calc ul a ti ons 

4. 1. 1 ARIEL 1 MAIN BODV 

Although the data in Figure 9 is "bumpy," the 0 = 84® profile for the main body 
indicates a nUUor central peak or bump, of height about 80 percent of ambient, 
within the depleted wake region where the minimum fs about 50 percent of ambient. 

Figure 10 shows transverse profiles computed for the Wake Of tiie main body. 

Of n^ (normalized ion density), n^ (normalized electron density), and^ (dimension- 
less potential) in the wake region downstream. The parameter values in Figure 10 
are = -4, M = 4. Xjj = l/lO. Thirteen major iterations (PoiSSon- Vlasov cycles) 
were computed. The profiles are in transverse planes at various distances down- 
stream. and all lengths are normalized by the body radius. Thus, z denotes axial 
distance downstream, in radii, with Z = 0 defined as the rear surface of the pillbox 
(looking into ttie wake); and r denotes radial or transverse distance from the axis 
(r = 1 is the body radius). The profiles Of n^. n^, and ^ are arranged vertically in. 
ondei^ Of increasing axial distance z. There are 8 values of z. namely, z = 0. 2. 

0. 5, 1. 6, 2. 0, 8. 6, 4. 0, 5. 0, and 6. 0. Each profile is constructed of 9 values 
Of r, namely, t » 0, 0. I, 0.3. 0.6, 0.8, 1.0, 1.2, 1.5. amd2.0. The Outer- 
boundary conditions are applied at Z = 8 and r 2; fon the main-body problem, the 
boundary condition at z « 8 is the fixed one, While floating conditions are used 
elsewhere. The profiles consist Of Straight-line segments connecting the values Of 
the functions (nj, n^, or ♦) computed at 12 grid points in the wake region. 

The features of the wake structure are as follows. The near wake (z < 1) is 
clearly depleted of both ions and electrons, with the ivn density lower than the 
electron density. Further downstream the wake becomes increasingly filled ih, 
between about z » l and z = 4, where z = 4 is the ton-Mach-number of radii 
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Figul'e 10. Ariel 1 Main Wody Wake Prdllles. ^ = -4, M = 4, 

Xjj » (10)’» = 0. 1 

downstream. In this range of z the iOn profiles tend to be relatively nolSy, Indica 
ting sensitivity to nximerleal errors, which may in turn imply a tendency toward 
physleal InStSbllUy. 'ftie trend of the loii and electron profiles suggests a radlally- 
Inward bulk motion of the plasma within the wake as if it were a fluid wave propaga 
ting Inward, piling up near the center, and bouncing out again as It moves down- 
stream. The disturbance has essentially died away at a * fi radii downstream. 



There Is ho well’-deflhed eehtral butnp aiiiiilair td that In the eJcpefimental 
data in Figure 0. 

4 , 1.2 ARIEL 1 ION PROBE 

Accbrdlhg to HehderSon artd Samir the profile in Figure 9 at fl = fiO® samples 
the ion probe wake etruoture. This structure is similar to that of the main body, 
having a below -ambient certtral peak within a depleted region of about the same 
width as that associated with the main body. The iOn probe iS about 1 . 1 Debye 
lengths-irt radius and is biased at about 28 kT (7 volts) negative with respect to 
space. Henoe it may be ejtpected to produce at least a pronounced focusing effect, 
as is borne out by the following computed results. 

Figure 11 shows transverse profiles computed for the wake of the ion probe. 

The notation is the same as that of Figure 10. The parameter values are = -28, 

M = 4, and.^Q = 1/1.7. Ten major iterations were computed. The outer boundary 
at 2 = 10 was placed Sufficiently far downstream to ensure that the disturbances of 
interest are contained within the grid. Moreover, a floating condition is employed.... 
there as well as elsewhere. (The number of grid points was larger..th.an in the 
main-body problem. ) 

The main features are as follows. The ion profiles at z = 1 and z = 2 show 
that a strong focusing effect occurs near the body. Further downstream, however, 
the disturbance dies Out; there is essentially nOne at z = 8 and beyond. 

The radially- inward and outward bulk motion Of the plasma as it fills the wake 
is again a fluid-like feature. Again, there is ilo persistent peak at the center of 
the wake as indicated by the data. If this were, an isolated body, the Henderson- 
Samir data would imply that the peak persists far downstream to beyond a distance 
of 33 radii. The present theoretical calculation indicates no Structure at z = 8 and 
beyond. 

4. 1. 3 COMMENT ON COMPARISON OF THEORY WITH EXPERIMENT 

A central bump may perhaps be generally expected on the basis of nonrigorous 
theorJBtiCal arguments invoking (1) electrostatic focusing effects Or (2) convergence 
of ion streams during the filling of the wake principally at a Mach number of radii 
downstream, or (3) a combination of these. However, previous theoretical 
calculations indicating sUch bumps have been deficient in some respect with regard 
to their rigorous applicability (for example, cold ions, infinitely- long cylinders, 
non-self-cortsistent). Slmilai. ly, there have been laboratory-simulation experiments 
\«rhich have Indicated bumps. However, it is pi'esShtly still di^hcult to simulate 
loft transverse velocity distributions IM-the laboi*atory, arid the effective ion 
temperature is generally too low. Herice, thsre does not exist thus far an unambigu- 
ous explariation of the Ariel 1 data. (Note; Tj was not measured, but was assumed 
here to be cquril to . Theo^etically, bumps in the wake have been predicted for 
Tj less than 
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It ts also of interest to note that a la^ge body produces a central ‘-core erihariee- 
merit, as will be showrt later, 

4/2 tixplorer 31 SutellUe 

The results of this secttori were obtained in the process of computlhg a Uutnber 
of solutions to be compared with in-Situ data obtained on the Explorer 31. The 
parameters adopted were suggested by Samir (private communication) based on 
8. different passes of the Explorer 31 satellite, as listed in his paper with Jew. 

A Small portion of two of these cases will be discussed here, without a quantitative 
Comparison with data, in order to illustrate a specific point. 

Fl^re 12 shows computed electron and ion density transverse profiles in the 
very near wake of a body with the parameters 

Xp = (6. 9)-l 
M =3.4 , 

■These are the parameters In the ease of Curve No. 1 of Samir and Jew. Tiie 
computed profiles in- Figure 12 are at z - 0. 2 radii downstream, that is, Similar 
to the lowest profiles in Figure 10 where the parameters are Of similar order. 

Here the vertical Scale (normalized density) is logarithmic, as opposed to Figure 
10. where it is linear. The ion densities are denoted by circles and the electron 
densities by squares. The principal features Shown in Figure l2 are as follows. 
For r greater than about 0. 8 radii, the ion density is higher than the electron 
density. Moreover, the ion density drops more abruptly in the vicinity of the 
"shoulder" (r^ 1) than the electron density. In the central wake both densities are 
far below normal, with the lOn density about an order of magnitude below the 
electron density. This is the usually- expected picture of near-wake structure. 

Figure 13, on the other hand, shows corresponding computational results for. 
the case of Curve No. 4 of Samir and Jew, where the parameters are: 



^0 “ -S '* 

Xj3 = (3. 1)-* 

M = 1. 1 

H(gre, the modt Significant change in the lower Mach number <1. 1). The curves in 
Figure 13 are tiualltatiVely different from thdhe id Figure 12. Thb new electron 
dat;sity profile id shifted downward slightly, but the hew ion density profile ts 
mdv6d up tb approicimately a constant roughly dttual to the nortnal valub. 
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The fact fcal the Ion «ake denalty la aho/e the electron aenetty at all r may be 
u«np.ctea from the point ot vie. ol "traditional" .ake theory, but aeema reaao.mble 
on the baaia ot Langmuir probe theory According m probe 
negative eloetroatatlc probe In a plaama ..U have adiacent to It a aheath In .h eh 
the ion denalty ekoeade the electron denalty. It the probe begtna to move alowly 
relattve to the plaama. One ettpecld the sbeath atructure at tlrat to be only aligbay 
Cbanged. «lth a continuation ot the prodomtnaoc# ot the Iona over the electrona. 

At auttlolently large velocity, however, the traditional wake atruWore .lib electron 
domination over Iona ahould appear aa In Figure 12. The value ot ton Mach number 
at «hlch the traoeltlon ahould occur haa not been predicted but can be eatabliahed 
by additional computations of the present type. 


5. \ IiARriE-BOIlY PROBLEM 

In this section vie consider the wake of a large body. 100 Debye lengths in 
radius. ^ The body IS in the form of a disk oriented normal to the flow. For this 
case (Figures. 14 and 15) the parameter values are 


Figure 14 



M = 4 . 


Figure 15 

= .4 

Xjj = (lOO)"^ 
M = 8 . 


Here the parameter values differ qualitatively from those of the preceding problems 
in that is So small. This Size Of moving body is larger than has been previous y 
treated by trajectory-following, that is. realistic, calculations. The reSuUs show 
what may be expected for the wake structure of large bodies in general. This case 
requires more effort (computer time and judicious selection of numerical para- 
meters) than that of a smaller body. 'Hie solutions Shown, therefore, are intended 
to be illustrative rather than accurate. 

SIX iterations, or PoisSOn-VlaSOv cycles. Were computed using the lon-de si y 
optlbb. lb «blch ducccdSlVe Iterated were bol mixed, dtartlbg .Itb Ibe beulral lob 
dendlly as ao Iblllal gueee. Tbe nominal number ot trajectories. 512. was ueed at 

all grid points. . j 

The proltlea ot n,. n.. and d in Figure 14 are cbnStrucled In Ihe same .ay and 

at the same grid points as In Figure 10. Tabulaled values are given In Parker. 

nie *ake la essentially "empty" of both tons and electrons bobvecn a - 0 and a 1. 

and begins to ttU up between a ■ 2 and a ■ 3. In this way. the ..ke la gual.lallvely 

similar to that in Figure 10. 
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Two Bfeti of Idn-donslty profiles ai*e shown on thfc left side of Flgu**e 14, the 
unlabeled prr riles for the fith Order ffith iteration), and the profiles labeled "A" 
for the 5th order. Comparison of the n^ -profiles with the 5th order Oj -profiles 
(labeled "A" to denote that the d-profUes and n^ -profiles in the figure are derived 
from these) indicates that the quasineutrality assumption is valid everywhere out- 
side a cone-shaped region rear the wake surface; the cone height along the axis is 
between 1 and 2 radii. This is in accord with expectation for a large body. Near 
the wake surface, however, quastneutrality is violated because the effective Debye 
length is large. The similarity of the n^-profiles (labeled "A") and the n^ profiles 
in Figure 14 is a consequence of near-quaSineutrallty. 

Comparison Of the 5th and 6th order ni-profiles (labeled "A” and unlabeled) 

in Figure 14 show that the solution IS reasonably converged for x = 1 and below, 

but that there is incomplete convergence at z = 2 and beyond. The Incomplete 
convergence and apparent structure at z = 2 and beyond may be artifactual due to 
insufficient numerical accuracy, (No attempt was made to achieve high accuracy 
since this was regarded as a preliminary run, ) The structure and lack of conver- 
gence arc seen to extend past z = 5, so that the downstream boundary should be 
placed further than at z = 6. 

Despite possible inaccuracies, one may infer additional physical conclusions 
indicated by Figure 14, namely, (l) the suggestion of a core of high (approximately 
ambient) density of ions and electrons on the axis, and (2) the occurrence at a 
potential well in the near wake, defined as a region with d-values below -4. The 
shading in ttte two lowest 0-profileS denote cross -sections of this well. The wake- 
surface normalized fluxes are 1. 1 XlO’® (5th order) and 2. 4 X lo"”^ (6th order) 
for ions, and 4. 3 X lo’® tr - electrons. The electron current density iS less tian 
eXp(-4), as would be expected In the presence of a potential well. 

The region of wake disturbance probably extends more than 6 radii downstream, 
and between 2 and 3 radii in the transverse direction. 

Another large-body case (Figure 15) is Similar to the previous large-body 
case except that the Mach number is increased from M = 4 to M = 8. Ten itera- 
tions were computed in which successive iterates were used without mixing, 
starting with uniform ambient ion density. (The latter starting condition was in- 
advertently different from that of the M = 4 calculation which was started with the 
neutral loh density, but this difference should become unimportant after many 
iterations. ) Similar statements may be made about the incompleteness of the 
convergence as in the M = 4 case. The 9th and lOth order ion densities are 
labeled "A " and unlabeied, respectively. On comparing these, the convergence 
seems fairly good at Z = 0. 5 and z = 1. Again, the disturbance extends beyond 
Z = 5, so that the downstream boundary should be moved further tnan z = 6. 


Despite tnBccuraclos. the conslstenoy is such that physical conclusions may 
be drawn as follows. In this case, the wake IS seen to felnaln empty further down 
stream than in the M M case. In addition, the suggostion la much stronger that 
there is a central core of ambient density for both lens and electrons along the 
axis. Moreover, the potential well is wider and longer than in thcM n 4 ease, 
although the depth is about ihe same. The normalized Wake-surface fluxes are 
7. 4 X 10- (9th order) and 4. 2 X io-^O (lOth order) for lorts. and 3. 7 X lo’3 for 

electrons. The electron flux Is slightly less than the M = 4 value, and is again 
less than exp (- 4 ). 

The conical region behind the disk where guaslneutrality breaks down Is now 
longer than in the M = 4 case, extending to between z = 4 and z = 5 along the axis. 
The region of wake disturbance Is probably longer than fi radii downstream. 

as in the M = 4 case, but may not extend beyond about 2 radii In the transverse 
direction. 
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